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ABSTRACT 

We modeled the V-band light curve of /3 Lyr with two stellar components plus an optically 
thick accretion disc around the gainer assuming a semidetached configuration. We present 
the results of this calculation, giving physical parameters for the stars and the disc, along 
with general system dimensions. We discuss the evolutionary stage of the system finding the 
best match with one of the evolutionary models of Van Rensbergen et al. According to this 
model, the system is found at age 2.30 x 10^ years, in the phase of rapid mass transfer, the 
second one in the life of this binary, in a Case-B mass-exchange stage with M = 1.58 x 10^^ 
Mq yr~^. This result, along with the reported rate of orbital period change and observational 
evidence of mass loss, suggests that the mass transfer in /? Lyrae is quasi-conservative. The 
best model indicates that (3 Lyrae finished a relatively large mass loss episode 31400 years 
ago. The light curve model that best fit the observations has inclination angle i = 81°, Mi = 
13.2 M0, Ma = 3.0 Mq, Ri = 6.0 R© and R2 = 15.2 R©. The disc contributes 22% to the 
total F-band light curve at quadrature, has a radius of 28.3 R© and the outer edge thickness is 
1 1.2 R©. The light curve model is significantly better with two bright regions in the disc rim 
with temperatures 10% and 20% higher than the disc outer edge temperature. We compare 
our results with earlier studies of this interacting binary. 

Key words: stars: early-type, stars: evolution, stars: mass-loss, stars: emission-line, stars: 
variables-others 



1 INTRODUCTION 

P Lyrae (Sheliak) is one of the best studied interacting binaries and 
also one of the more enigmatic bright objects in the sky. The liter- 
ature for this object is extensive and only brief highlights will be 
given in this short introduction. The interested reader is directed 
to the reviews by Sahade (1980) and Harmanec (2002). The cur- 
rent picture of the system includes a B6-B8II donor star transfer- 
ring mass onto an early B-type gainer. The hotter star is hidden 
by an optically and geometrically thick disc (Huang 1963; Wilson 
1974; Hubeny & Plavec 1991; Skulskii 1992, Linnell 2000) and 
some peculiar features can be explained by the presence of a bipo- 
lar jet emanating from the disc (Harmanec et al. 1996; Hoffman et 
al. 1998, Ak et al. 2007). Interferometry of /3 Lyrae showed for the 
first time the direct image of a gravitationally distorted star (Zhao 
et al. 2008; Schmitt et al. 2009). Non-conservative evolution was 
inferred from the presence of a radio nebula surrounding the binary 
system (Umana et al. 2000). The system has an orbital period of 
12.9 days increasing at a rate of 19 s yr~^ (Harmanec & Scholz 
1993). Among the notable properties of this system are: (1) the un- 
explained non-orbital 282.4 d photometric periodicity, (2) the rever- 
sal of the depth of primary and secondary eclipses and their gradual 
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evanescence towards the ultraviolet and (3) the presence of six sys- 
tems of infrared, optical and ultraviolet lines characterized by their 
morphology and radial velocity behavior (e.g. Harmanec 2002). 



Many solutions for the system parameters have been proposed 
in the literature, some of them including an accretion disc around 
the gainer. In this paper we offer a complementary view of the sys- 
tem, and apply our light curve (LC) synthesis code to fit the pub- 
lished l/-band light curve, obtaining a new set of representative 
system parameters. After obtaining these parameters, we search for 
a self-consistent evolutionary model by exploring a grid of theoret- 
ical evolutionary tracks for binary stars. The synthetic models are 
those by Van Rensbergen et al. (2008, 2011); they consider brief 
episodes of liberal or non-conservative evolution. The paper is or- 
ganized as follows. In Section 2 we describe the light curve synthe- 
sis code including an accretion disc around the gainer and present 
the results of the light curve analysis. In Section 3 we explore the 
evolutionary stage of f5 Lyr obtaining age and mass transfer rate 
and a theoretical set of system parameters. In Section 4 we provide 
a discussion of our results and a comparison with previous investi- 
gations. In Section 5 we give our conclusions. 
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Table 1. Parameters of the best fit model for the /3 Lyrae V-band fight- 
curve obtained by solving the inverse problem considering an accretion disc 
around a gainer in critical rotation regime. 



Quantity Quantity 





2852 


M \MrA 
' h l-^ ' W J 


13.16 ± 0.3 


S(0 — C)'^ 




A/f \ A/f^l 


9 q7 -1- n 9 




0.0415 




6.0 ± 0.2 


'I J 


86 1 + 5 




15.2 ± 0.2 




n QA + n OS 


iogioi/h 


4.0 ±0.1 


Td[K] 


8200 ± 400 


logl09c 


2.5 ±0.1 


de [dorb] 


0.192 ± 0.009 


bol 


-6.3 ±0.2 


do [dorb] 


0.01 ±0.01 


bol 


-4.7 ±0.1 


ax 


3.8 ±0.3 


lorbPo] 


58.5 ±0.3 


fh 


20.2 ±0.5 


7^d[R0] 


28.3 ±0.3 


Fh 


1.00 


de[R0] 


11.2 ±0.2 


Th[K] 


30000 ± 2000 


dc[R0] 


0.6 ±0.1 




1.21 ±0.1 


logioTh 


4.48 


^hsn 


16.2 ±2.0 


logioTo 


4.12 


Ahs[°] 


324.6 ± 7.0 


logioLh 


4.42 




-13.5 ±5.0 


logioLc 


3.81 


Abs = Tbs/T(j 


1.12 ±0.1 






SbsT] 


53.8 ±5.0 






AbB[°] 


107.3 ±9.0 






Ah,c,d 


0.10, 0.35, 0.59 








12.12 ± 0.02 






Ho 


2.296 ±0.02 







Fixed parameters: q = Mc/Mh = 0.226 - mass ratio of the com- 
ponents, Tc = 13300K - temperature of the less-massive (cooler) donor, 

Fc = 1.0 - filling factor for the critical Roche lobe of the donor, Fjj = 
Rh/Rzc = 1.0 - filling factor for the critical non-synchronous lobe of 
the hotter, more-massive gainer (ratio of the stellar polar radius to the crit- 
ical non-synchronous lobe radius along z-axis for a star in critical rotation 
regime), fc = 1.00 - non-synchronous rotation coefficients of the donor, 
/3h c = 0.25 - gravity-darkening coefficients of the components. 

Note: n - number of observations, S(0 — C)^ - final sum of squares of 
residuals between observed (LCO) and synthetic (LCC) light -curves, fJrms 
- root-mean-square of the residuals, i - orbit inclination (in arc degrees), 
F(j = R(j/Ryc - disc dimension factor (the ratio of the disc radius to the 
critical Roche lobe radius along y-axis), T^j - disc-edge temperature, de, 
dc, - disc thicknesses (at the edge and at the center of the disc, respectively) 
in the units of the distance between the components, ax - disc tempera- 
ture distribution coefficient, /h - non-synchronous rotation coefficient of 
the more massive gainer (in the critical rotation regime), - tempera- 
ture of the gainer, A^s.bs = bs/Td - hot and bright spots' temper- 
ature coefficients, 0hs.bs and A^s bs " spots' angular dimensions (radius) 
and longitudes (in arc degrees), S^ad " angle between the line perpendicular 
to the local disc edge surface and the direction of the hot-spot maximum 
radiation, Ajj c d " albedo coefficients of the system components and the 
accretion disc, fib ^ - dimensionless surface potentials of the hotter gainer 
and cooler donor, [A^q], TZ^ ^ [Rq] - stellar masses and mean radii 
of stars in solar utiits, logioffh.c " logarithm (base 10) of the system com- 
ponents effective gravity, Af^^j - absolute stellar bolometric magnitudes, 
<^orb [R©]; ^d[R©], dcfR©], dcfR©] - Orbital semi-major axis, disc ra- 
dius and disc thicknesses at its edge and center, respectively, given in solar 
units, logioTh^c, logioLh,c " logarithm (base 10) of the system compo- 
nents effective temperature and luminosity. 

2 LIGHT CURVE ANALYSIS 

2.1 Model for an optically thick disc around the gainer 

Here we give a brief description of the disc model that we ap- 
ply to /3 Lyrae. The basic elements of the binary system model 
with a plane-parallel disc and the corresponding Ught-curve synthe- 



sis procedure are described in detail by Djurasevic (1992a, 1996). 
The model and code have been widely used, slightly improved and 
tested during our recent research of intermediate-mass interacting 
binaries (e.g. Djurasevic et al. 2010, 2011, 2012, Garrido et al. 
2012, Mennickent et al. 2012). 

We assume that the disc in /3 Lyr is optically and geometri- 
cally thick. The disc edge is approximated by a cylindrical surface. 
In the current version of the model (Djurasevic et al. 2010), the 
vertical thickness of the disc can change linearly with radial dis- 
tance, allowing the disc to take a conical shape (convex, concave 
or plane-parallel), i.e. the disc presents a triangular cross-section 
at both sides of the star. The geometrical properties of the disc are 
determined by its radius (Rd), its vertical thickness at the edge (de) 
and the vertical thickness at the center or inner boundary (dc). 

The cylindrical edge of the disc is characterized by its 
temperature, T^, and the conical surface of the disc by a radial 
temperature profile obtained by modifying the temperature distri- 
bution proposed by Zola (1991): 

Tir)=T, + in-Ta)[l-{-^)r- (D 

We assume that the disc is in physical and thermal contact with 
the gainer, so the inner radius and temperature of the disc are equal 
to the temperature and radius of the star (Rh, Th). The temperature 
of the disc at the edge (T4) and the temperature exponent (ar ), 
as well as the radii of the star {Rh) and of the disc {Rd) are free 
parameters, determined by solving the inverse problem (see Section 
2.3). 

The model of the system is refined by introducing active re- 
gions on the edge of the disc. The existence of these active regions 
is supported by hydrodynamical simulations of /3 Lyrae as we will 
discuss in Section 4.1. The active regions have higher local tem- 
peratures so their inclusion results in a non-uniform distribution of 
radiation. The model includes two such active regions: a hot spot 
(hs) and a bright spot (bs). These regions are characterized by their 
temperatures Ths,bs angular dimensions (radius) 9hs,bs and longi- 
tudes \hs,bs- The longitude A is measured clockwise (as viewed 
from the direction of the -i-Z-axis, which is orthogonal to the orbital 
plane) with respect to the line connecting the star centers (-i-X-axis), 
in the range 0-360 degrees. These parameters are also determined 
by solving the inverse problem. 

Due to the infall of an intensive gas-stream, the disc surface 
in the region of the hot spot becomes deformed as the material ac- 
cumulates at the point of impact, producing a local deviation of ra- 
diation from the uniform azimuthal distribution. In the model, this 
deviation is described by the angle 6rad between the line perpen- 
dicular to the local disc edge surface and the direction of the hot 
spot maximum radiation in the orbital plane. 

The second spot in the model, i.e. the bright spot, approxi- 
mates the spiral structure of an accretion disc, predicted by hydro- 
dynamical calculations (Heemskerk 1994). The tidal force exerted 
by the donor star causes a spiral shock, producing one or two ex- 
tended spiral arms in the outer part of the disc. The bright spot can 
also be interpreted as a region where the disc significantly deviates 
from the circular shape. For further information about these active 
regions in the model, see Djurasevic et al. (2012). 

One limitation of the code at its present implementation is the 
lack of a detailed treatment for the donor irradiation by the disc 
part facing the donor, including hot spot. This effect could be po- 
tentially important in close binaries with a large difference in stel- 
lar and disc temperatures. However, it should be a second-order ef- 
fect compared with the stellar and disc flux contributions (including 
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donor irradiation by tlie gainer) already implemented in our code. 
This is demonstrated in the very good fit obtained with the orbital 
light curve. 

Another limitation is the lack of consideration of a "third 
light", a source of continuum emission that is not eclipsed; in the 
case of (5 Lyr it could be the jet reported by Harmanec et al. (1996). 
However, studies show that this third light has a minor influence in 
the total flux of P Lyrae, being larger in the ultraviolet, but much 
less evident in the optical (e.g. Linnell 2000). Hence, it is reason- 
able to assume that the optical flux is dominated by the stars and 
the disc. 

2.2 The V-hsmA time series 

We have used the V-band light curve and non-linear ephemeris of 
P Lyrae published by Harmanec et al. (1996). This datasct consists 
of 2582 V-band magnitudes obtained by differential photometry in 
a 36 year interval at 12 different stations. In most stations the re- 
ported mean rms error per one observation of the corresponding 
check star is mostly between 0.01 and 0.03 mag. The interested 
reader can find more details including the procedure of homoge- 
nization of the multi-station light curve in Harmanec et al. (1996). 

Since the long photometric cycle of 282.4 d has a low amph- 
tude of 0.024 mag, comparable to the error of an individual data- 
point, viz. 0.027, and this long cycle is probably not strictly periodic 
(Harmanec et al. 1996), we decided to work with the original V- 
band light curve without intending to remove the long cycle, which 
could introduce additional uncertainties. 

2.3 The fitting procedure 

The light-curve fitting was performed using the inverse-problem 
solving method based on the simplex algorithm, and the model of 
a binary system with a disc described in the previous section. We 
used the Nelder-Mead simplex algorithm (see e.g. Press et al. 1992) 
with optimizations described by Dennis and Torczon (1991). While 
the direct problem comprises the calculation of the light curve from 
model parameters given a priori, the inverse problem is the process 
of finding the set of parameters that will optimally fit the synthetic 
light curve to the observations. For more detail see e.g. Djurasevic 
(1992b). 

To obtain reliable estimates of the system parameters, a good 
practice is to restrict the number of free parameters by fixing some 
of them to values obtained from independent sources. Thus we 
fixed the spectroscopic mass ratio to g = 0.226 and the donor tem- 
perature to T2 = 13300 K, based on published values (e.g. Har- 
manec 2002, Linnell 2000). Our own LC simulations made with 
the parameters of the best solution but changing the mass ratio 
confirmed this q value. In addition, we set the gravity darkening 
coefficient and the albedo of the gainer and the donor to I5h,c = 
0.25 and Ah,c = 1.0 in accordance with von Zeipel's law for radia- 
tive shells and complete re-radiation (Von Zeipel 1924). The limb- 
darkening for the components was calculated in the way described 
by Djurasevic et al. (2010). 

The possible values of free parameters are constrained by im- 
posing the lowest and highest values which seems reasonable based 
on previous studies of this binary. Here are the ranges for the fitted 
parameters: 

• Inclination: 70 to 90 degrees. 

• Disk dimension factor (the ratio of the disk radius and the ra- 
dius of the critical Roche lobe along the y-axis): 0.5 to 1.0. 



• Disk edge temperature: 6000 to 12000 K. 

• Disk thickness: 0.01 to 0.3 distance between the components 

Oorb- 

• The exponent of the disk temperature distribution: 1 to 5. 

• Non-synchronous rotation coefficient for the gainer (in critical 

rotation regime): 5 to 25. 

After the first fit, these ranges are decreased according to the 
results of the first iteration. 

We treated the rotation of the donor as synchronous, i.e. we 
used the non-synchronous rotation coefficient, defined as the ratio 
between the actual and the Keplerian angular velocity, fc = 1.0. 
Since it is assumed that the donor has filled its Roche lobe (i.e. 
the filling factor of the donor was set to Fc = 1.0), the synchro- 
nization is expected to be the consequence of the tidal effects. In 
the case of the gainer, however, the accreted material from the disc 
is expected to transfer enough angular momentum to increase the 
rate of the gainer up to the critical velocity as soon as even a small 
fraction of the mass has been transferred (Packet 1981, de Mink, 
Pols & Glebbeek 2007). This means that the gainer fills its corre- 
sponding non-synchronous Roche lobe for the star rotating in the 
critical regime, and its dimensions and the amount of rotational dis- 
tortion are uniquely determined by the factor of non-synchronous 
rotation. For (3 Lyr we assumed critical rotation for the gainer, and 
estimated a non-synchronous rotation factor fh = 20.2 in the crit- 
ical rotation regime. Our additional modeling, obtained consider- 
ing the synchronous rotation regime for the gainer, does not differ 
much from the critical rotation case. 

We were able to model the asymmetry of the light curve very 
precisely by incorporating two regions of enhanced radiation on the 
disc: the hot spot (hs), and the bright spot (bs). The hot spot and the 
bright spot in our model are located on the edge side of the disk and 
are described by the longitude of the center of the spot, the angular 
radius of the spot, and the temperature ratio of the spot and the 
unperturbed local temperature of the disk. The difference between 
the temperature of the spot and the local unperturbed temperature 
of the disk is what results in the difference in brightness. The hot 
spot can be located at longitudes between 270 and 360 degrees, 
and the bright spot can be located at any longitude. The angular 
radius of the spots was constrained to the range from to 30 degrees 
for the hot spot, and from to 70 degrees for the bright spot; the 
temperature ratio for the spots can be from 1.0 to 1.5. 

2.4 Results of the light curve fitting 

Compared to the model without the active regions on the disk, the 
model which includes the hot and bright spots (hs+bs) fits the ob- 
servations much better (see Fig. 1). The dominant feature of our 
optimal model is the bright spot, and the parameters describing its 
size, temperature and location are the most important for improving 
the fit. The existence of a bright spot indicates significant deviation 
of the disk from cylindrical shape, or other phenomena that are not 
directly modeled, such as a spiral arm in the outer part of the disk. 
The parameters of the hot spot don't have as much influence on the 
quality of the fit because the bright spot is significantly larger in 
size. 

An attempt to model the system without including the active 
regions on the edge of the accretion disk showed that the observed 
light curve cannot be described well with such a model. This is 
especially noticeable in the secondary minimum and around the 
secondary maximum where there is a significant asymmetry in the 
light curve (see Fig. 1). The 0-C residuals between the observed 
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Figure 1. Observed (LCO) and model's synthetic light-curves (without spots - no sp, with hot spot - hs, with bright spot - bs and with hot spot and bright 
spot - hs+hs on the accretion disk edge) of ,6 Lyrac obtained by analyzing photometric observations and the corresponding final 0-C residuals between the 
observed and model's synthetic light curves. The O-C plots show also the final S(0 — C)^ for the corresponding model and polynomial fit of the residuals 
presented by the solid line. In the right-hand column a thin horizontal line shows magnitude. 



("LCO") and synthetic light curve ("no sp") clearly show the defi- 
ciencies of this model for optimal fitting of the observations. 

If a hot spot region is added to the model of the disk, the model 
fits the observations better in the secondary minimum, but still does 
not describe the asymmetry in the region of the secondary maxi- 
mvun, and introduces problems with the fitting of the primary mini- 
mum. The third model we tested has only a bright spot on the edge 
of the disk and gives a somewhat better fit, but with significant dif- 
ferences between the model and the observations in certain phase 
intervals (see Fig. 1). Finally, Fig. 1 shows the optimal model of the 



system, which contains both the hot spot and the bright spot regions 
(hs-i-bs) on the disk. The picture clearly shows that this model can 
successfully fit the observations and the light curve asymmetries 
that were problematic in other models we considered. 

Apart from this graphical presentation, the advantages of our 
optimal model can be demonstrated by applying the Fisher sta- 
tistical F-test, which shows that the models without spots (S'l = 
S(0 - Cf = 5.338) and the models with a single (hot - & = 
E(0 - Cf = 5.313 or bright - S3 = S(0 - Cf = 5.244) spot 
can be discarded with the probability of 95% when compared to the 
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Figure 2. The dependence of S(0 — C) '^ on basic model parameters. Each panel shows the variation of one parameter while all others are icept fixed on the 
optimal values obtained for the best-fitting model, which includes the hot and bright spots (hs+bs) on the accretion disk. The dashed line shows the minimum 
value of S(0 - C)2 = 4.904 in each panel. 
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optimal model (hot spot and bright spot - Sopt ~ S(0 — C) = 
4.904). Namely, the mentioned values give the following quanti- 
ties: Fi = Si/Sopt = 1.0885, F2 = S2/Sopt = 1.0834 and 
F3 = S3 /Sopt = 1.0693, while the Fisher statistical F-test for the 
probability of 95% gives the critical value of the F-distribution F= 
1.064. This clearly shows that the model with active hot and bright 
regions on the disk fits the observations much better than models 
without spots (Fl) and with single hot spot (F2) or bright spot (F3) 
because these ^"1,2,3 calculated from the data are greater than the 
critical value of the F-distribution for probability 95%, as can also 
be seen in Fig. 1. On the other hand, adding more bright spot re- 
gions does not improve the quality of the fit. 

Fig. 2 shows how varying of basic model parameters around 
the optimal values (obtained for the best model, which includes 
the active regions on the disk) affects the E(0 — C)^ value. Each 
panel shows the variation of one model parameter while all others 
are kept fixed on their optimal values. The scale in all plots is the 
same, to make easier to note which parameters influence the quality 
of the fit at most. 

The model is relatively insensitive to the temperature of the 
gainer and to the exponent of the disk temperature distribution. This 
is to be expected since the disk is seen edge-on, and it covers a large 
part of the gainer, so the inner regions of the disk don't contribute 
much to its total radiation. The value found by us, viz. = 30.000 
K, is consistent with that used by other authors (e.g. Hubeny & 
Plavec 1991). The model also shows relatively low sensitivity to 
changes in disk radius (parameter Fd), inner disk thickness (dc), 
and the non-synchronous rotation coefficient of the gainer (fh)- It 
is more sensitive to changed in the orbital inclination, the outer disk 
thickness (de) and temperature (T^). 



2.5 The best light-curve model 

The fit, O — C residuals, individual donor, disc and gainer flux 
contributions and the view of the optimal model at orbital phases 
0.25, 0.50 and 0.75, obtained with the parameters estimated by the 
light curve analysis, are illustrated in Fig. 3. We note that residuals 
show no dependence on orbital or long-cycle phases, except a larger 
random scatter around main eclipse, a feature already evident in 
the light curve. We find that the best fit model of f3 Lyrae contains 
an optically and geometrically thick disc around the hotter, more 
massive gainer star (Table 1). With a radius of Rd « 28.3Rq, 
the disc is 4.7 times larger than the central star (Rh w 6.QRq). 
The disc has a convex shape, with central thickness rfc ~ O.6-R0 
and edge thickness de ~ 11.2Rq. The temperature of the disc 
increases from Ta = 8200 K at its edge, to Th = 30000 K at the 
inner radius, where it is in thermal and physical contact with the 
gainer. The relatively large disc temperature gradient explains the 
big difference between disc thickness at the inner and outer edges. 

In the best model the hot spot with 16° angular radius covers 
9% of the disc outer rim and it is situated at longitude \hs ~ 325°, 
roughly between the components of the system, at the place where 
the gas stream falls onto the disc (Lubow & Shu 1975). The tem- 
perature of the hot spot is approximately 21% higher than the disc 
edge temperature, i.e. Ths ~ 9900 K. Although including the hot 
spot region into the model improves the fit, it cannot explain the 
light- curve asymmetry completely. By introducing one additional 
bright spot, larger than the hot spot and located on the disc edge 
at \bs ~ 107°, the fit becomes much better. This bright spot has 
Tbs « 9200 K and with 54° angular radius covers 30% of the disc 
outer rim. 

A Fourier analysis (Lenz & Breger 2005) and a Phase Disper- 
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Figure 3. Observed (LCO) and synthetic (LCC) light curves of /3 Lyr, ob- 
tained by analyzing photometric observations; final 0-C residuals between 
the observed and optimum synthetic light curves; fluxes of the donor, the 
gainer and of the disc, normaUzed to the donor flux at phase 0.25; the views 
of the optimal model at orbital phases 0.25, 0.50 and 0.75, obtained with 
parameters estimated from Ught curve analysis. We indicate the mass ra- 
tio (q), the approximate positions of the inner Lagrangian point (Li), the 
hot spot (hs), the bright spot (6s), and the gainer non-synchronous rotation 
factor defined in the text (fh)- 
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sion Minimization analysis (Stellingwerf 1978) of residuals of the 
best fit model confirm the long periodicity of 282 days found by 
Hannanec et al. (1996) with the same dataset but with other light 
curve model. 



3 ON THE EVOLUTIONARY STAGE OF S LYR 

In this section we compare the system parameters obtained in 
Section 2.2 and listed in Table 1 with those predicted by binary 
evolution models including epochs of non-conservative evolution. 
We inspected the 561 conservative and non-conservative evolu- 
tionary tracks by Van Rensbergen et al. (2008) available at the 
Center de Donnees Stellaires (CDS) looking for the best match for 
the system parameters found for p Lyr. Models with strong and 
weak tidal interaction were studied, although only the latter ones 
should allow critical rotation of the gainer. Following Mennickent 
et al. (2012), a multi-parametric fit was made with the synthetic 
{Sij,k) and observed (Ok) stellar parameters mass, temperature, 
luminosity and radii (obtained from the light-curve fitting), and 
the orbital period, where i (from 1 to 561) indicates the synthetic 
model, j the time tj and k (from 1 to 9) the stellar or orbital 
parameter. Non-adjusted parameters were mass loss rate, Roche 
lobe radii, chemical composition, fraction of accreted mass lost by 
the system and age. For every synthetic model i we calculated the 
quantity Xi,j at every tj defined by: 

Xlj = il/N)EkWk[{Si,j,k - Ok)/Ok]^ (2) 

where N is the number of observations (9) and Wk the statistical 
weight of the parameter Ok, calculated as: 

Wk = y/Ok/e{Ok) (3) 

where e{Ok) is the error associated to the observable Ok- The 
model with the minimum corresponds to the model with the best 
evolutionary history of (3 Lyr. The absolute minimum Xmin identi- 
fies the age of the system along with the theoretical stellar and or- 
bital parameters. The high accuracy of the orbital period dominates 
the search for the best solution in a single evolutionary track (Fig. 
5), but the others parameters play a role when comparing tracks 
corresponding to different initial stellar masses. 

We find the absolute minimum (viz. 0.0291) in the weak 
interaction model with initial masses of 12 and 7.2 Mq and ini- 
tial orbital period of 2.5 days; their parameters are shown in Table 
2. All other models give larger x^ and were rejected as reliable 
solutions. A comparison of x^ for these models with the optimal 
solution indicates that the best fit stands out among other solutions 
in the multi-parametric space (Fig. 6). The solution very near to the 
absolute minima in Fig. 6 (easily visible in the bottom panels) cor- 
responds to the strong interaction model for the same initial condi- 
tions; it has x^ = 0.0308 and their parameters, also given in Table 
2 for reference, differ by a small and almost always neglectable 
amount from the adopted parameters. 

The corresponding evolutionary tracks for the gainer and 
donor stars are shown in Fig. 4, along with the position for the 
best model for /S Lyr. We observe a relatively good match for the 
gainer star parameters but a small mismatch with the donor tem- 
perature and luminosity. The observed donor is slightly fainter and 
cooler than the model prediction. The difference amounts to 0.9% 
in log T and 2.7% in log L. This small discrepancy remains when 
considering the next more reliable solution, those with strong tidal 



interaction and whose parameters are given in Table 2. We pro- 
pose that this mismatch could be due to a limitation of the existing 
grid of synthetic models available, that are not dense enough in ini- 
tial orbital period and initial masses to give account of all possible 
evolutionary combinations, and eventually could not reproduce the 
13 Lyrae evolution perfectly. For instance, only initial masses for 
the donor of 7, 8, 9, 12 and 15 M© are considered in the range 
compatible with j3 Lyrae. Whereas this discretization of initial pa- 
rameters probably does not produce a large discrepancy with the 
actual system parameters (as suggested by the absence of compe- 
tent solutions near the absolute minima in Fig. 6), it is likely that 
small differences between observed and synthetic parameters, as 
observed with the donor, occur. At present, we cannot quantify this 
effect. 

The best fit indicates that /3 Lyr is found inside a burst of mass 
transfer, the second one in the hfe of this binary (Fig. 5). The donor 
is an inflated (J?c = 14.4 Rq) and evolved 2.2 M© star with its core 
completely exhausted of hydrogen. According to the best model the 
system has now an age of 2.3 X lO'^ yr and M = 1.58 x 10"^ M© 
yr~^. The system reverted its mass ratio not at the present semide- 
tached stage, but about 8 Myr ago, during a much stronger episode 
of mass transfer. The small difference between Mh and in Ta- 
ble 2 is a purely numeric effect due to the iterative process of cal- 
culation; the transferred mass equals the accreted mass in the long- 
term, i.e. according to the best model P Lyrae is in a conservative 
mass transfer stage. An inspection of the best-model evolutionary 
track indicates that: (i) the system is found 147.000 years after 2nd 
Roche lobe contact and (ii) in this interaction episode the system 
has lost (previously) 0.03 M© into the interstellar medium and the 
gainer has eventually accreted 1.55 M©. The last figures were de- 
rived from the mass balance among two epochs in the evolutionary 
track: the epoch corresponding to the start of mass loss and the 
present epoch. However, if critical rotation impedes the accretion 
of more mass (something not considered by the binary evolution 
codes), part of this mass could have been ejected into the interstel- 
lar medium or remains still stored in the accretion disc. 



4 DISCUSSION 

In this section we compare our /3 Lyrae parameters with those pre- 
viously reported in the literature. We also give a semblance of the 
future appearance of the system, just after finishing the mass trans- 
fer process, based on the prediction of the best synthetic model. 



4.1 Comparison with published system and disc parameters 

A range of stellar masses, inclinations and mass ratios have been 
published for (3 Lyrae during the past 50 years, but the latest calcu- 
lations seems to indicate gainer and donor masses around 13 M© 
and 3 M© , respectively (Harmanec 2002). Our stellar masses, listed 
in Table 1, are consistent with these numbers. 

Similarly to other authors (e.g. Hubeny & Plavec 1991, here- 
after HP, and Linnell 2000), we have found that a geometrically 
thick disc significantly improves the fit to the light curve in /3 Lyrae. 
In our model the disc contributes 22% to the total F-band light 
curve at quadrature, has a radius of 28.3 R© (HP find 25.5 R© and 
Linnell 30 R©) and the thickness goes from 0.58 R© at the inner 
edge to 1 1.2 R© at the rim (HP find 0.34 R© to 14 R©, respectively, 
and Linnell finds 16 R© at the rim). 

Our derived mass transfer rate M = 1.58 x 10~® M© yr~^ 
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Figure 4. Evolutionary tracks for the binary star model from Rensbergen et al. (2008) that best fit the data. Donor (right black track) and gainer (left gray 
track) evolutionary paths are shown, along with the parameters derived from the LC fit. The best fit is reached at the time corresponding to the model attached 

to the axis by dashed lines, that is characterized in Table 2. The mismatch for the donor is discussed in the text. Stellar sizes are proportional to the circle 
diameters. An expanded view is shown in the right panel where half of the points are shown for the donor track for best visualization. 



matches the canonical value commonly used in numerical simula- 
tions (1-4 X 10"® M0 yr~\ e.g. Bisikalo et. al. 2000, Nazarenko 
& Glaztmova 2003). This canonical value is often calculated by 
assuming conservative mass transfer and the observed value for 
the orbital period change (19 s yr~^; Harmanec & Scholz 1993). 
In fact, a binary interchanging material in a conservative way has a 
constant period change given by: 



(Huang 1963). For /3 Lyr this equation gives M = (2.2 ± 0.2) x 
10~® Mq yr~^, where the major contribution to the error budget 
comes from the donor mass uncertainty quoted in Table 1. 

We notice two important issues here: (1) The best fit model 
is conservative, but observations (discussed in the next paragraph) 
suggest mass loss at the present epoch and (2) our M value is 
near the conservative figure derived from Eq. 4. Whereas the first 
point could indicate that the models still do not reproduce the ac- 
tual physics of mass loss processes, the second point suggests that 
the mass loss rate probably is much smaller than the mass transfer 
rate, producing a mass transfer practically conservative or, in other 
words, quasi-conservative. 

Jets and nebular radio emission have been considered signa- 
tures of mass loss in /3 Lyrae. From the evolutionary track discussed 



in Section 3, we have estimated that 0.03 Mq have been lost in the 
present episode of Roche lobe overflow. Umana et al. (2000) es- 
timate 0.015 Me from the analysis of the radio nebula but using 
M ~ 10^^ Mq yr^^ for the donor wind ionized by the radia- 
tion field of the hotter companion. An earlier estimate based on 
the mass exchange models by De Greve and Liimell (1994, here- 
after DGL94) is 0.43 Mq. Our mass loss estimate is linked to the 
nature of the mass loss process considered in the Van Rensbergen 
et al. (2008, 201 1) models. These models consider mass loss driven 
by radiation pressure from a hot spot located in the stream impact 
region on the stellar surface or accretion disc edge. The mass loss 
extracts angular momentum from the system, specifically from the 
gainer. The best model indicates that /3 Lyrae already passed by a 
non-conservative epoch 3 1400 years ago, and now it is in a conser- 
vative regime. The fact that observations actually show some mass 
loss suggests a quasi-conservative stage instead, as stated previ- 
ously. 

We find two bright regions along the disc rim at longitudes Xhs 
= 325°and Xbs = 107°; the first one 1.21 times hotter than the disc 
rim, with 16° angular radius covering 9% of the outer disc rim, and 
the second one 1.12 times hotter than the disc rim, with 54°angular 
radius, covering 30% of the outer disc rim. These regions might 
be tentatively identified with shock regions revealed in hydrody- 
namical simulations of mass transfer in P Lyrae by Nazarenko & 



On the accretion disc and evolutionary stage of p Lyrae 9 



14 
12 
10 

Is 

^6 

o 

^ 4 
2 


10 



1 r 



0,1 r 



0,01 



-1 

>. 

1-2 h 

f -3 
o 

5,-4 - 

O 

1 -6 



I I I I I I I I I I I I I I I I I I I I I I I I I I l|l I I I I I I I I N 



I I I I I I I I I I I I I I I I I I I I I I I I I I I 



I I I I I I I I I I 



12 14 16 18 20 22 24 26 
age (Myr) 

Figure 5. Mc, (defined in Eq. 2) and Mc for the best evolutionary track. 
The vertical dashed line indicates the position for the best model. 



Glazunova (2006) and in close binaries in general by Bisikalo et al. 
(1998, 1999, 2003). 

Interestingly, the hot spot has also been inferred from spec- 
tropolarimetric observations by Lomax et al. (2012). They report 
a hot spot at the disc outer edge with a projected linear extension 
between 22 and 33 Rq, equivalent to angular sizes 27.4 and 60.8 
degrees (7.6 and 16.9 % of the outer disc rim), depending on the 
photometric bandpass. The center of the hot spot is roughly at or- 
bital phase 0.47, near the position of the hot spot found by us at 
orbital phase 0.40 {Xhs = 325°). 



4.2 Comparison with previous evolutionary studies 

We find /3 Lyrae in a Case-B mass transfer stage with an age of 2.30 
X 10^ yr. The donor has exhausted hydrogen in its core and the 
gainer has a hydrogen concentration X^h = 0.43. The parameters 
for the best non-conservative evolutionary model are given in Table 
2. The late evolutionary stage of the donor is consistent with the 
spectroscopic finding of He enrichment and extreme CNO cycling 
by Balachandran et al. (1986). 

The evolutionary stage of fi Lyrae was previously investigated 
by DGL94. These authors determined that a conservative solu- 
tion led to several discrepancies with observed data. They find that 
P Lyrae is reaching the end of the mass transfer, having evolved in 
a non-conservative way from an initial system of 9 M© -i- 7.65 M© 
and orbital period 4 days (our best model starts with 12 M© -i- 7.2 



M0 and orbital period 2.5 days). Their model parameters are given 
in Table 2 along with those of our best model. 

In general, the stellar parameters derived from our evolution- 
ary solution compares relatively well with those of DGL94. How- 
ever, significant differences can be observed in the initial masses 
and the initial orbital period. In addition, DGL94 model gives a 9 
Myr older binary, a slightly more inflated donor and a mass transfer 
rate 75% larger than ours. 

The above differences are not negligible and point to some 
physical difference between both models. One of the possible rea- 
sons for the differences might be the more massive binary consid- 
ered in our solution, along with the shorter initial orbital period. 
This means a faster evolution (a younger system as observed) but 
also that non-conservative effects should be stronger than in the 
case studied by DGL94. Another possible reason is the way both 
models deal with mass loss, usually parametrized with the param- 
eter /3, defined as the mass fraction accreted by the gainer, so 1-/3 
is the mass fraction lost by the system. DGL94 assume /3 constant 
through the whole mass transfer process whereas Rensbergen et al. 
consider P variable according to the physics of their model. In any 
case, we caution that we cannot go further interpreting the differ- 
ences between both models since the mass loss processes are still 
not well understood, neither included with all physical details in the 
models. At this stage, these differences can be interpreted as upper 
estimates for the accuracy of the initial binary parameters, its age 
andM. 



4.3 The future of /? Lyrae 

According to the best synthetic model presented in Section 3, the 

donor will lose contact with its Roche lobe after 72100 yr, becom- 
ing a detached binary, with stellar components of 1.78 M© and 
13.77 Mq. The orbital period will be 23.47 days and the stars will 
be in basically the same present evolutionary stage. The model does 
not cover the posterior evolution of the gainer, but if the disc and 
rapid gainer rotation still persist at that time, the system might be 
easily classified like a binary Be star. This conclusion depends on 
how the angular momentum and circumstellar mass are distributed 
during the interaction epochs, something that future models need 
to consider carefully. 



5 CONCLUSIONS 

In this paper we fit the published optical light curve of j3 Lyrae with 
a light curve synthesis code based on the simplex algorithm, con- 
sidering stellar plus disc components. In addition, the evolutionary 
stage of the binary was studied by comparing our best system 
parameters with those of grids of synthetic evolutionary tracks. 
Hence, we have made use of the published data in an entirely new 
and complementary way. A comparison with earlier studies was 
also given. The main results of our research are: 

• The V light curve can be modeled with an intermediate mass 
semidetached binary including an accretion disc around the gainer. 
The best system parameters including stellar and disc temperatures, 
masses, luminosities and sizes are given in Table 1. 

• Our study of the evolutionary stage of d Lyrae indicates that it 
is in a Case-B mass transfer stage, with age 2.30 x 10^ years and 
M = 1.58 X 10~^ M0 yr~^. The donor has exhausted hydrogen 
in its core and the gainer is slightly evolved with central hydrogen 
fraction Xch = 0.43. 
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Table 2. The parameters of the Van Rensbergen et al. (2008, labeled vR08) models that best fit the f} Lyr data of Table 1 along with those of the model by De 
Greve & Linnell (1994, labeled DGL94). Subindexes "s" and "w" refer to the strong and weak interaction models for the same set of initial parameters, where 
the weak case is the best solution followed by the strong case, as discussed in the text. The hydrogen and helium core mass fractions (Xc and Yc) are given 
for the cool and hot star. Errors are given when available, corresponding approximately to half of the grid step at a given parameter. 



quantity 


vROSu, 


vR08s 


DGL94 


quantity 




vR08s 


DGL94 


age (yr) 


2.30070(4)E7 


2.30457(1)E7 


3.1924E7 


period (d) 


12.95(9) 


12.95(8) 


13.23 


(M0) 


2.248(6) 


2.231(6) 


2.95 


Mh (Me) 


13.308(6) 


13.251(6) 


13.27 


Mc (Mq yr-i ) 


-1.58(1 l)E-5 


-1.38(2)E-5 


-2.77E-05 


Mh (Mq yr-i) 


1.44(2)E-5 


1.47(4)E-5 


2.57E-5 


log Tc (K) 


4.162(1) 


4.161(1) 


4.135 


log Tu ( K) 


4.482(1) 


4.482(1) 


4.444 


log Lc (L0) 


3.917(1) 


3.911(1) 


3.917 


log Lf, (Lq) 


4.442(1) 


4.449(1) 


4.329 


Rc (R©) 


14.36(5) 


14.29(3) 


16.26 


Rh (R© ) 


6.023(2) 


6.076(2) 


6.33 


^CC 


0.00(0) 


0.00(0) 


0.000 




0.43(0) 


0.42(0) 


0.479 




0.98(0) 


0.98(0) 


0.980 




0.55(0) 


0.56(0) 


0.501 




• The best evolutionary model indicates that the system is found 
147.000 years after 2nd Roche lobe contact. In this interaction 
episode the system has lost 0.03 M© into the interstellar medium 
and the gainer might have accreted 1.55 M©. A fraction of this 
mass could still remain in the disc. 

• We provide arguments supporting the view that /3 Lyrae is in a 
quasi-conservative evolutionary stage, having finished a relatively 
large mass loss episode 31400 years ago. By quasi-conservative. 



we mean a stage with a systemic mass loss much less that the mass 

transfer rate. 

• Discretization of the grid of synthetic models could explain 
the small discrepancy observed in donor temperature. At present 
we cannot quantify this effect, but argue it should be small. 

• We find significant evidence for two bright regions in the 
disc rim with temperatures 10 and 20% higher than the disc outer 
edge temperature. They are located at longitudes A = 107° and 
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325°respectively. The cooler bright spot region covers 30% of the 
disc outer rim and the hotter one 9%. These hot regions might be 
tentatively identified with those shock regions revealed by hydro- 
dynamical simulations of the mass transfer process. 
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